Statistical approach of the modulational 
instability of the discrete self-trapping 

equation 

Anca Vi§inescu, D. Grecu 

Department of Theoretical Physics 
National Institute for Physics and Nuclear Engineering 
"Horia Hulubei" 
P.O.Box MG-6, Magurele, Bucharest, Romania 
e-mail: avisin@theorl . theory, nipne. ro 
dgrecu@theorl . theory, nipne. ro 

Abstract 

The discrete self-trapping equation (DST) represents an useful 
model for several properties of one-dimensional nonlinear molecular 
crystals. The modulational instability of DST equation is discussed 
from a statistical point of view, considering the oscillator amplitude 
as a random variable. A kinetic equation for the two-point correla- 
tion function is written down, and its linear stability is studied. Both 
a Gaussian and a Lorentzian form for the initial unperturbed wave 
spectrum are discussed. Comparison with the continuum limit (NLS 
equation) is done. 
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1 Introduction 



The discrete self-trapping (DST) equation 

dCLf2 , \ I 1 2 / \ 

i— - u a n + \{a n+1 + a n -i) + fi\a n \ a n = (1) 

is a typical equation for a system of harmonically coupled nonlinear oscilla- 
tions Q, || relevant for several physical problems. We mention here only 
Davydov's model of energy transport in a-helix structures in proteins [Q], 
||, where (0) appears as a certain approximation of the model. In ([!]) a n is 
the complex classical dimensionless amplitude of the oscillator of frequency 
ujq in the n-th molecule, and A, fi (of dimension of frequency) are the coupling 
constants between nearest neighbour oscillators and the one-site nonlinearity 
respectively. It is well known that depending upon of the parameters and 
the chosen initial condition the equation (j^) can lead either to self-trapping 
(i.e. local modes or solitons), or to chaos, or to a mixture of the above two 
behaviours [[jj, f2[, [|]. Instead of (]l|) we shall consider the equation 

l ~dT + ^ an+1 + a ' a ~ 1 ' + ^l ™' a " = ' > 

which is obtained if a n —>■ a n e~ tu)at . This equation admits plane wave solutions 
with constant amplitude 

(the lattice constant is taken equal with unity) but with an amplitude de- 
pending dispersion relation 

ui{k) = — 2A cos k — n\a\ 2 

This is a Stokes wave solution and it is well known to be unstable at small 
modulation of the amplitude (Benjamin- Feir or modulational instability) ||- 
|| . The aim of this note is to study the modulational instability of equation 
(0) from a statistical point of view, considering (Xyi clS db random variable. In 
doing this we shall follow the procedure used by several authors to discuss the 
effects of randomness on the stability of weakly nonlinear waves, especially 
in hydrodynamics P|, RIO . 
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In the next section a kinetic equation for a two-point correlation function 
will be obtained. Using a Wigner-Moyal transform the equation is written 
in a mixed configuration-wave vector space. The linear stability around a 
homogeneous basic solution is discussed in section 3. An integral stability 
equation is derived, very similar with the dispersion relation of the linearized 
Vlasov equation in ionized plasmas. Two forms for the spectrum of the 
initial unperturbed condition will be considered, namely a Gaussian and 
a Lorentzian form and in the limit of vanishingly small bands widths the 
increment of the modulational instability is calculated. Comparison with 
the continuum limit, when ([!]) transforms into the nonlinear Schrodinger 
equation is done. Few concluding remarks are also presented. 

2 Kinetic equation for two-point correlation 
function 



Introducing the displacement operator by a n ±\ = ea™a n the equation (Q) 
becomes 

da d 

i— - + 2Acosh— a n + n\a n \ 2 a n = 0. (3) 
at on 

In order to find a kinetic equation we write (§) for n — m, multiply it by 
add it to the complex conjugated of @ for n = n 2 multiplied by a ni and 



a 



a-, 1 



finally take an ensemble average. One obtains 

^ < «m< > +2 A (cosh - cosh < a ni a* n2 > + 

which beside the two point correlation function p(ni, n 2 , t) =< a ni {t)a* n2 {t) > 
contains also four-point correlation functions. If a n corresponds to a Gaus- 
sian process, and this property is retained during the evolution, a four-point 
correlation function factorizes exactly in products of two-point correlation 
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< a ni a* ni a ni a* n2 >= 2 < a ni a* 2 >< a ni a* m >= 2p(n 1 ,n 2 )a 2 (n 1 ) (5) 

where a 2 (n) =< a n a* n > is the ensemble average of the mean square am- 
plitude. Although the factorization (Kf) is true only for a Gaussian process 
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we shall assume to be at least approximately valid also for processes slightly 
different from a Gaussian one, and it represents the main approximation of 
the present analysis. 

It is convenient to use a Wigner-Moyal transform |TI| . One introduce the 
new variables 

M— , m = n 1 -n 2 . (6) 

Then the equation @ becomes 

% % + 4A sinh \m sinh ^ + Ma2{M + y } " a2{M " y ))p = (7) 

We consider a chain of N molecules and impose cyclic boundary condi- 
tions. The Fourier transform of the two-point correlation function is defined 
by 

777 777 

F(k, M,t) = J2 e- tkm p(M + -, M - -, t) (8) 

m 11 

where k takes values in the first Brillouin zone (BZ), k G (— 7T, n). The inverse 
formula is 

mm 1 P^. 1 r n 

p(M+- , M — —,t) = -J2 e lkm F(k, M,t) = — e* km F(k, M, t)dk. (9) 

For m = one obtains 

a\N, t) = ^ E F(k, M,t) = ±- f W F(k, M, t)dk. (10) 

1\ k Z7T J-tt 

Now Fourier transforming equation (|7|) we get 

dF ^ • , • , 1 9 r, 

— — h 4A sm /c smh — — — t + 
dt 2 DM 

which is the expected nonlinear evolution equation for F(k, M, t) in a mixed 
configuration- wave number space (M, k). Using the definition (§) we see that 
F(/c, M, t) is a periodic function in the reciprocal space, F{k + 2n) = F{k). 
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3 Stability analysis 



As the unperturbed problem we shall consider a basic solution Fo(k) inde- 
pendent of M and t. This is the random counterpart of the Stokes wave in 
a deterministic approach. A small perturbation around this homogeneous 
background is considered, namely 

F{k,M,t) = F {k) + eF 1 (k,M 7 t) (12) 



According to ([10]) we have also 

a 2 (M, t) = a 2 + eof (M, t) (13) 

where 

al = ^-T F (k)dk 

af(M,t) = — f F 1 (k,M,t)dk (14) 

When (^) is introduced into ([TT]) , neglecting terms of order e 2 , the following 
linear evolution equation for i*\ is obtained 

^1 ^ , i 1 d n 

_ + 4AsmJfe S mh-^ ? F+ 

- (-iy +1 d^F d^aj(M) 
i^x (2j - 1)!2 2 J'- 1 dP^ 1 dA^- 1 1 J 

Looking for a plane wave solution 

F x (k,M,t) = h{k)e i{KM - nt) 

after little algebra the following stability integral equation is found 

* 4A sm ~ 

The modulational instability is related to Q complex with a positive imagi- 
nary part, ImQ > 0. It is convenient to compare with the similar result 
for the continuum case of the nonlinear Schrodinger equation || 

t-gfe (17) 
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Although there are significant differences between the two expressions, when 
the width of the spectrum Fo(k) is vanishingly small, the final results will 
look very similar. 



3a. Gaussian spectrum 

As a first example let us assume F Q (k) to be a Gaussian function 



2ti 



k 



F (k) = a\e~^. (18) 

a 

This expression doesn't satisfy the periodicity condition but for a vanishingly 
small the errors introduced are negligible. Also the relation ( |T4| ) is satisfied 
up to exponentially small terms. 

It is convenient to introduce the new integration variable t = ~^{k ± y) 
and the notations 

i / n k\ 

z± = —^[— y -±tan— (19) 

V2aV2Asinir 2 J v ' 

1 / K 
f±(t) = —j=- ( sin \[2ot ± tan — (1 — cos Vzcrt) 
V2o" V 2 



Then ( fL6|) becomes 

" e~* 2 f — f ^—r j tft = 1. (20) 



2 /. ^ 



27rt A sin If ./--§- ^2+ - /+ 2_ — /_ 

In leading order in a the integral (p0|) can be evaluated using the steepest 
descent method [|T3[| . Denoting Cr±(i) = t 2 + Zn(z± — f±(t)), t± the zeros 
of the first derivatives = 0, A± = \ ^2t f° r t = t±, and extending 

the integration limits to infinity the integral is given by 

7T ( J= e -G + (t + ) _ _j_ e -G-{t-)\ . (21) 



In the limit a -C 1 we have approximatively 

v 2Asin A 2 

Then the integral becomes 



t± ^ t- = *[% a * e - G± ('±) ~ ^ and A± ~ 1. 



-2 v / 27rcxtanf 

2Ashi^) 2 ~ ( taI1 f ) 2 
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Considering f2 purely imaginary, f2 = zOj, we finally get 



a = 4Asinyya^-sin 2 y (22) 

and an instability is obtained (fij > 0) if and A > and if sin 2 | < a§ t. 
3b. Lorentzian spectrum 

A simpler example is a Lorentzian form for -Fq(^) 



Fo(k) = al Pv ^f 2 . (23) 
sm ^ + 



It satisfies the periodicity condition and relation (|14|). The unperturbed 
two-point correlation function is easily calculated using (|2l|) in the definition 
relation (|T§D. Straightforward calculations give 
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Po(m) = irTWTTFTrtF (24) 

representing an exponentially decreasing law. For p < 1 we have po{m) ~ 
a 2 e- 2pm . 

In order to calculate the integral ([16]) it is convenient to introduce the 
new integration variable t = tan y. Then the integral is over the whole real 
axis and can be done in the t-complex plane. In the new variable Fo(k ± y) 
writes 

F ( k± *L)^ a 2 2pVTTJ?(l + t 2 ) 

01 2 } (1 + cosf + 2p 2 )t 2 ±2(sinf)t+(l-cosf + 2p 2 ) 1 } 



having poles at 
where 



t\ 2 = —a ±ib t 12 = a±ib 



K_ 

2 

1 + cos f + 2p 2 ' 1 + cos f + 2p 2 



sin f 2p^l+p 2 



(26) 
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Considering f2 purely imaginary, f2 = iVLi and denoting z = ., Qi k we have 

4A sin 

also 

1 1 + t 2 

— z- 



, , . -r — sin k zt 2 + 2it + z 

4A sin ~ 

having poles at 



t 3 = % 14 = —1 

z z 

We shall consider z as a small quantity and consequently £4 3> 1. Closing the 
contour in the lower complex half-plane i its contribution can be neglected 
in the first order. Therefore we shall take into account only the poles 
and after straightforward calculations the relation (|T6|) becomes 



fia 2 MA + MX 
Asinf X 2 + M 2 



(27) 



where 



A = 1 + a 2 - b 2 , B = 2ab 

X = zA + 2b, M = 2a — zB. (28) 



When p< 1 we approximate 



sin tt 2p 



a - , , 'jr . & - T-T^Tk (29) 



1 + cos y 1 + cos -j 



terms of order p 2 being neglected. Then ( p7|) can be considerably simplified 
and finally give us 



,\ K I V-? ■ 2 K o 1 + cos T + cosk 1 
^ = 4A sin -\\l^a 2 - sin 2 — - 2p - g ) . ( . J >( I j 

Modulational instability occurs for A and fj, > 0, sin 2 y < an d if V is 



smaller than a critical value. This result is similar with a previous one |L6 
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obtained with a simplified Lorentzian form for Fo(k). Both results (B^) and 
(|30D can be compared with the similar results obtained in the NLS case ||. 



(G) 



^ ) = 2tfa*|W-^a§-^--p| (31) 




where the superscript G/L refers to Gaussian/Lorentzian form of Fo(k). It 
is easily seen that ([31]) are obtained in a long wave limit (K <C 1). In the 
Lorentzian case both relations (^) and ( |31~1) show a behaviour similar with 
the well known phenomena of Landau damping in plasma physics [15 



namely with increasing of p the imaginary part f2j can become negative and 
no instability develops. 

In conclusion a complete discrete discussion of the randomness effects on 
the MI of the self trapping equation was done and the discrete effects are 
easily seen in the final results, compared with the similar ones found for the 
NLS equation. 

Helpful discussions with Dr. A.S. Cdrstea are kindly acknowledged. This re- 
search was supported under the contract 66, CERES Program, with the Ministry 
of Education and Research. 
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